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Abstract 


New  models  of  non-stationary  processes  have  been  introduced.  We  showed 
that  the  covariance  of  locally  stationary  processes  can  be  modeled  as  pseudo- 
differential  operators.  The  spectrum  of  such  processes  is  estimated  by  ap- 
proxiinating  the  Karhunen-Loeve  basis  with  a  local  cosine  basis,  that  is  op- 
Umized  with  a  best  basis  search  algorithm.  Applications  to  geophysics  have 
been  studied.  Locally  dilated  processes  are  a  different  kind  of  non-stationary 
processes  that  appear  in  image  processing  and  in  physical  phenomena  that 
involve  Doppler  effects.  It  was  shown  that  the  dilation  parameters  of  such 
processes  can  be  estimated  with  a  wavelet  transform,  through  the  solution 
of  a  partial  differential  equation  in  the  scale-space  plane.  An  application 
concerns  the  reconstruction  of  three  dimensional  surfaces  from  texture  gra¬ 
dient  m  images.  The  last  part  of  this  project  was  devoted  to  the  analysis  of 
the  distortion-rate  function  of  wavelet  image  transform  codes.  By  modeling 
imaps  as  element  of  Besov  spaces,  we  calculated  precise  analytical  formula 
of  distortion-rates,  which  are  verified  by  numerical  experiments. 


1  Introduction 

Most  processes  encountered  in  physics  or  signal  processing  are  non-stationary 
but  are  often  approximated  by  stationary  processes  because  of  the  lack 
o  models  and  estimation  procedures  to  analyze  non-stationarity.  In  the 
hrst  part  of  this  project  we  have  introduced  a  general  scheme  to  estimate 
the  covariance  of  non-stationary  processes  by  searching  for  an  approximate 
Karhunen-Loeve  basis  among  a  dictionary  of  bases.  This  general  scheme 
was  applied  to  locally  stationary  processes.  A  model  of  locally  station¬ 
ary  processes  based  on  pseudo-differential  operators  was  introduced  We 
showed  that  the  spectrum  of  these  locally  stationary  processes  can  be  effi¬ 
ciently  estimated  with  a  dictionary  of  local  cosine  bases  [1,2].  Applications 
to  geophysics  have  been  studied. 

Locally  dilated  processes  is  a  second  important  class  of  non-stationary 
processes.  These  processes  appear  in  image  textures,  in  finance  and  in  phys¬ 
ical  phenomena  including  Doppler  effects.  Such  processes  are  created  by 
dilating  or  compressing  a  stationary  process,  with  scaling  factors  that  vary 
wi  h  time  (or  space).  The  relevant  information  is  often  stored  in  these  scal¬ 
ing  factors  that  must  be  estimated.  We  proved  that  the  covariance  of  these 
dilated  processes  satisfy  a  transport  differential  equation  in  a  scale-space 
fficTore"  ^  ^  based  algorithm  was  derived  to  estimate  the  scaling 

Wavelet  image  transform  codes  are  currently  the  most  performant  for 
ill  image  compression,  yet  there  is  no  theoretical  understanding  of  their 
distortion  rate  properties  because  images  are  realizations  of  highly  non- 
s  ationary  processes  that  are  not  well  understood.  We  showed  that  the 
current  distortion-rate  theory  does  not  apply  to  understand  the  properties 
o  hese  codes.  We  introduced  a  new  approach  to  compute  the  distortion 
ra  e  unction,  by  estimating  the  decay  rate  of  sorted  wavelet  coefficients 
This  IS  equivalent  to  use  models  based  on  Besov  spaces. 


2  Covariance  Estimation  and  Locally  Stationary 
Processes 

Second  order  moments  characterize  entirely  Gaussian  processes  and  are  of- 
en  sufficient  to  analyze  stochastic  models,  even  though  the  processes  may 
t  be  Gaussian.  When  processes  are  wide-sense  stationary,  their  covariance 
defines  a  convolution  operator.  Many  spectral  estimation  algorithms  allow 
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one  to  estimate  the  covariance  operator  from  a  few  realizations,  because 
It  IS  diagonalized  with  Fourier  series  or  integrals.  When  processes  are  not 
sUtionary,  in  the  wide-sense,  covariance  operators  may  have  complicated 
time  varying  properties.  Their  estimation  is  much  more  delicate  since  we 
do  not  know  a  priori  how  to  diagonalize  them.  The  ideas  and  methods  of 
Calderon  and  Zygmund  in  harmonic  analysis  have  shown  that  although  we 
are  not  able  to  find  the  basis  which  diagonalizes  complicated  integral  op¬ 
erators  in  general,  it  is  nevertheless  possible  to  find  well  structured  bases 
which  compress  them.  This  means  that  the  operator  is  well  represented  by 
a  sparse  matrix  with  such  a  basis.  This  approach  allows  characterization 
of  large  classes  of  operators  by  the  family  of  bases  which  do  the  compres¬ 
sion.  We  have  shown  that  the  ability  to  represent  covariance  operators  by 
sparse  matrices  in  a  suitable  basis  leads  to  its  efficient  estimation  from  a  few 
realizations. 

In  collaboration  with  Professor  Papanicolaou  and  Dr.  Zhang  [1]  we  have 
concentrated  our  attention  on  the  class  of  locally  stationary  processes,  that 
is,  processes  whose  covariance  operators  are  approximately  convolutions. 
Since  cosines  and  sines  diagonalize  the  covariance  of  stationary  processes, 
it  is  natural  to  expect  that  local  cosine  functions  are  “almost”  eigenvectors 
of  locally  stationary  processes.  This  property  is  formalized  by  postulating 
that  the  covariance  operator  is  well  approximated  by  a  nearly  diagonal  one 
in  an  appropriate  local  cosine  basis.  We  have  shown  that  if  the  covariance 
operator  is  a  pseudo-differential  operator  of  a  specified  class,  then  the  process 
is  locally  stationary. 

To  estimate  the  covariance  operator  of  a  locally  stationary  process  we 
search  for  a  local  cosine  basis  which  compresses  it  and  estimate  its  matrix 
elements.  The  size  of  the  windows  of  a  suitable  local  cosine  basis  must 
be  adapted  to  the  size  of  the  intervals  where  the  process  is  approximately 
stationary.  Since  we  do  not  know  in  advance  the  size  of  approximate  sta- 
tionarity  intervals,  we  have  introduced  an  algorithm  that  searches  within  a 
class  of  bases  for  a  “best”  basis,  to  compress  the  covariance  operator.  This 
search  IS  done  using  data  provided  by  a  few  realizations  of  the  process.  It 
minimizes  the  Hilbert  Schmidt  norm  of  the  coefficients  along  the  diagonal 
of  the  covariance  matrix.  One  can  indeed  prove  that  the  Karhunen-Loeve 
basis  is  the  basis  that  minimizes  this  norm.  For  locally  stationary  processes, 
we  have  derived  a  fast  implementation  of  the  search  for  a  best  local  cosine 
basis  based  on  the  local  cosine  trees  of  Meyer,  Coifman  and  Wickerhauser. 

A  central  issue  was  to  prove  that  our  covariance  estimator  is  consistent 
when  the  sample  size  N  increases  to  infinity.  This  consistency  was  obtained 
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with  a  modified  algorithm  obtained  through  a  collaboration  with  Professor 
Donoho  and  Dr.  von  Sachs  [2].  The  central  idea  is  to  regularize  the  empirical 
covariance  coefficients  calculated  from  the  data  with  a  wavelet  thresholding 
algorithm.  We  proved  that  for  appropriate  locally  stationary  processes,  if 
the  size  of  the  maximum  eigenvalue  is  bounded  then  the  best  basis  covariance 
estimator  converges  to  the  true  covariance. 


3  Locally  Dilated  Processes 

Locally  dilated  processes  are  obtained  with  a  smooth  change  of  variable  7(#) 
of  a  stationary  process,  with  possibly  an  amplitude  modulation  a{t) 

The  derivative  Y{t)  can  be  interpreted  as  a  time  varying  scaling  factor.  Such 
processes  appear  in  image  processing  and  number  of  physical  problems  that 
incorporated  Doppler  effects.  In  image  processing,  a  homogeneous  texture 
that  can  be  modeled  with  a  stationary  process  X  will  appear  as  locally 
dilated  in  the  image  because  of  the  perspectivity  effect.  This  is  called  a 
“texture  gradient”.  Recovering  this  local  dilation  allows  one  to  compute  the 
t  ree  dirnensional  coordinates  of  the  surfaces  visualized  in  the  scene,  up  to 
a  multiplica,tion  factor.  In  collaboration  with  Maureen  Clerc  [3],  we  studied 
the  estimation  of  a{t)  and  j{t)  from  a  single  realization  of  y(t). 

We  suppose  that  the  variations  of  a(t)  and  j{t)  happens  on  a  larger  scale 
than  the  fine  scale  variations  of  X{t).  We  define  the  covariance 

Ciu,v)  =  E{Yiu+^)Y{u-^)}. 

The  scale  separation  hypothesis  allows  one  to  prove  that  C{u,v)  nearly 
satisfies  a  differential  equation 

^logC(u,t;)  _  dlogQ;(M)  dlogy(u)  d\ogC{u,v) 

du  du  d\o^v  ' 

This  equation  can  be  interpreted  as  a  transport  equation  in  scales  where 
(log7'(u)  j  is  a  velocity  along  scales.  We  show  that  we  get  nearly  optimal 

statistical  estimators  of  and  from  a  single  realization  of 

y  (t)  by  computing  the  covariance  in  a  wavelet  basis  at  fine  scales  [3].  The 
wavelet  transform  of  a  signal  /  is  defined  by 

=  f*°° 

J — oo  yd  d 
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The  variance  of  the  wavelet  transform  of  X{t)  is  R{b,a)  =  E{WY(b,a)]. 

At  fine  scale,  we  proved  that  the  wavelet  transform  satisfies  the  a  transport 
equation  similar  to  (1) 

d log R{u,  a)  ^  d\oga{u)  dlogY{u)  dlogR(u,a) 
du  du  5  log  a  ■ 

An  Iterative  algorithm  has  been  designed  to  estimate  7(t)  by  solving  (2)  from 
the  empirical  estimates  of  the  wavelet  transform  covariance  R(u,  a)  [3],  This 

algorithm  has  been  implemented  and  tested  numerically  on  one-dim"ensional 
signals. 


4  Understanding  Image  Transform  Codes 

Research  on  image  compression  led  to  a  culture  shock  between  several  math¬ 
ematicians  who  got  interested  in  this  topic  through  the  discovery  of  wavelet 
orthonormal  bases,  and  signal  processing  engineers  that  had  developed  effi¬ 
cient  algorithms  for  many  years.  If  one  can  approximate  an  image  with  few 
wavelet  image  coefficients,  many  mathematicians  thought  that  these  bases 
should  be  remarkably  well  adapted  to  perform  bit  compression.  Most  signal 
processors  have  replied  that  the  problem  is  much  more  difficult  since  it  also 
involves  a  careful  analysis  of  the  performance  of  scalar  quantizations  and 
bit  allocations.  If  the  signals  to  be  encoded  are  realizations  of  a  Gaussian 
process,  under  the  high  resolution  quantization  hypothesis,  we  know  nearly 
everything  on  the  performance  of  a  transform  code.  For  an  average  of  .R  bits 
per  pixel,  the  mean-square  error  D{R)  decays  proportionally  to  2“^*^  with 
a  constant  that  depends  upon  the  bit  allocation  and  the  basis.  The  optimal 
basis  IS  then  the  Karhunen-Loeve  basis,  and  wavelet  bases  are  a  priori  not 
pod  approximations  of  Karhunen-Loeve  bases  for  images.  It  is  known  that 
images  are  not  realizations  of  Gaussian  processes,  but  it  is  often  believed 
that  the  use  of  Gaussian  and  high  resolution  approximations  give  an  appro¬ 
priate  estimate  of  the  distortion-rate  curve  D{R).  Hence  the  widespread  use 
of  the  coding  “gam”  formula  which  is  only  justified  in  this  framework. 

Current  image  transform  codes  operate  below  1  bit  per  pixel.  For  such 
high  compression  rates,  we  have  shown  that  the  classical  transform  code 
theory  yields  a  wrong  estimate  of  the  distortion-rate  D{R)  [4,5].  If  one  opti¬ 
mizes  the  splar  quantization  and  bit  allocation  then  D{R)  depends  mostly 
on  ^e  precision  of  signal  approximations  with  few  non-zero  decomposition 
coefficients  in  the  orthogonal  ba^is.  Yet,  as  emphasized  by  signal  processors. 
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there  are  important  constants  involved,  that  depend  upon  the  encoding  ef¬ 
ficiency  of  the  positions  of  zero  coefficients.  At  low  bit  rates,  for  a  number 
of  bit  per  pixel  R  that  is  below  1,  we  proved  [4,5]  that 

D{R)  «  (1  -p  K)Do{j)  ~  (3) 

where  K  and  r  are  two  constants  that  are  specified,  and  Dq{z)  is  the  error 
when  approximating  a  signal  with  a  proportion  of  2  non-zero  basis  coef¬ 
ficients.  For  natural  images  decomposed  in  wavelet  or  cosine  bases,  the 
parameter  a  is  typically  of  the  order  of  1.  It  is  related  to  the  exponent  of 
the  smallest  Besov  space  that  includes  the  signal. 

Transform  code  algorithms  can  be  improved  by  an  embedding  strategy 
which  sends  first  the  larger  amplitude  coefficients  and  progressively  refine 
their  quantization.  The  improvement  provided  by  the  partial  ordering  of  em¬ 
bedded  codes  has  been  analyzed  mathematically  and  evaluated  numerically 
for  wavelet  zero-trees  [4,5]. 
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